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Abstract 

In this work, we have presented a simple analytical approximation scheme for 
generic non-linear FBSDEs. By treating the interested system as the linear decoupled 
FBSDE perturbed with non-linear generator and feedback terms, we have shown that 
it is possible to carry out a recursive approximation to an arbitrarily higher order, 
where the required calculations in each order are equivalent to those for standard 
European contingent claims. We have also applied the perturbative method to the 
PDE framework following the so-called Four Step Scheme. The method is found to 
render the original non-linear PDE into a series of standard parabolic linear PDEs. 
Due to the equivalence of the two approaches, it is also possible to derive approximate 
analytic solution for the non-linear PDE by applying the asymptotic expansion to 
the corresponding probabilistic model. Two simple examples are provided to demon- 
strate how the perturbation works and show its accuracy relative to known numerical 
techniques. The method presented in this paper may be useful for various important 
problems which have eluded analytical treatment so far. 
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1 Introduction 



In this paper, we propose a simple analytical approximation for backward stochastic dif- 
ferential equations (BSDEs). These equations were introduced by Bismut (1973) [JJ for 
the linear case and later by Pardoux and Peng (1990) [10] for the general case, and have 
earned strong academic interests since then. They are particularly relevant for the pricing 
of contingent claims in constrained or incomplete markets, and for the study of recursive 
utilities as presented by Duffie and Epstein (1992) [2]. For a recent comprehensive study 
with financial applications, one may consult Yong and Zhou (1999) [15], Ma and Yong 
(2000) [9] and references therein. 

The importance of BSDEs, or more specifically non-linear FBSDEs which have non- 
linear generators coupled with some state processes satisfying the forward SDEs, has risen 
greatly in recent years also among practitioners. The collapse of major financial institu- 
tions followed by the drastic reform of regulations make them well aware of the importance 
of counterparty risk management, credit value adjustments (CVA) in particular. Even in a 
very simple setup, if there exists asymmetry in the credit risk between the two parties, the 
relevant dynamics of portfolio value follows a non-linear FBSDE as clearly shown by Duffie 
and Huang (1996) [3]. We have recently found that the asymmetric treatment of collateral 
between the two parties also leads to a non-linear FBSDE [5j. Furthermore, in May 2010, 
regulators were forced to realize the importance of mutual interactions and feedback loops 
in the trading activities among financial firms, shocked by the astonishing flash crash of 
the Dow Jones index by almost 1, 000 points. Once we take the feedback effects from the 
behavior of major players into account, we naturally end up with complicated coupled 
FBSDEs. 

Unfortunately, however, an explicit solution for a FBSDE is only known for a simple 
linear example. In the last decade, several techniques have been introduced by researchers, 
but they tend to be quite complicated for practical applications. They either require 
one to solve non-linear PDEs, which are very difficult in general, or resort to quite time- 
consuming simulation. Although regression based Monte Carlo simulation has been rather 
popular among practitioners for the pricing of callable products, the appropriate choice 
of regressors and attaining numerical stability becomes a more subtle issue for a general 
FBSDE. In fact, in clear contrast to the pricing of callable products, one cannot tell if 
the price goes up or down when one improves the regressors, which makes it particularly 
difficult to select the appropriate basis functions. 

In this paper, we present a simple analytical approximation scheme for the non-linear 
FBSDEs coupled with generic Markovian state processes. We have perturbatively ex- 
panded the non-linear terms around the linearized FBSDE, where the expansion can be 
made recursively to an arbitrary higher order. In each order of approximation, the re- 
quired calculations are equivalent to those for standard European contingent claims. In 
order to carry out the perturbation scheme, we need to express the backward components 
explicitly in terms of the forward components for each order of approximation. For that 
purpose, we propose to use the asymptotic expansion of volatility for the forward com- 
ponents, which is now widely adopted to price various European contingent claims and 
compute optimal portfolios (See, for examples [Til HI H3 E] and references therein 
for the recent developments and review.). In the case when the underlying processes have 
known distributions, of course, we can directly proceed to a higher order approximation 
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without resorting to an asymptotic expansion. 

We have also studied a perturbation scheme in the PDE framework, or in the so- 
called Four Step Scheme [8], for the generic fully-coupled non-linear FBSDEs. We have 
shown that our perturbation method renders the original non-linear PDE into the series 
of classical linear parabolic PDEs, which are easy to handle with standard techniques. We 
then provided the corresponding probabilistic framework by using the equivalence between 
the two approaches. We have shown that, also in this case, the required calculations in a 
given order are equivalent to those for the classical European contingent claims. As a by- 
product, by applying the asymptotic expansion method to the corresponding probabilistic 
model, it was actually found possible to derive an analytic expression for the solution of 
the non-linear PDE up to a given order of perturbation. Therefore, our method can be 
interpreted as a practical implementation of the Four Step Scheme in the perturbative 
approach. 

The organization of of the paper is as follows: In Section 2, we will explain our new 
approximation scheme with perturbative expansion for generic decoupled FBSDEs. Then, 
in Section 3, we shall apply it to the two concrete examples to demonstrate how it works 
and test its numerical performance. One of them allows a direct numerical treatment by a 
simple PDE and hence it is easy to compare the two methods. In the second example, we 
will consider a slightly more complicated model. We compare our approximation result 
to the detailed numerical study recently carried out by Gobet et al. (2005) [6j using a 
regression-based Monte Carlo simulation. In Section 4, we explain how to use standard 
asymptotic expansion procedures to express the backward components explicitly when 
the forward components do not have known distributions. In Section 5, we will give an 
extension of our method to the fully coupled non-linear FBSDEs under the PDE frame- 
work, and then formulate the equivalent probabilistic approach in Section 5. Appendix 
contains slightly different scheme for coupled non-linear FBSDEs which may be useful for 
the actual application. 

2 Approximation Scheme 
2.1 Setup 

Let us briefly describe the basic setup. The probability space is taken as P) and 

T £ (0, oo) denotes some fixed time horizon. W t = (Wf, ■■■ , W[)* , < t < T is Revalued 
Brownian motion defined on ($7, J 1 ", P), and (-Pt){o<t<T} stands for P-augmented natural 
filtration generated by the Brownian motion. 

We consider the following forward-backward stochastic differential equation (FBSDE) 

dV t = -f(X t ,V t ,Z t )dt + Z t -dW t (2.1) 
V T = $(X T ) (2.2) 

where V takes the value in R, and Xf £ R rf is assumed to follow a generic Markovian 
forward SDE 

dX t = l0 (X t )dt + 7 (X t ) ■ dW t . (2.3) 

Here, we absorbed an explicit dependence on time to X by allowing some of its components 
can be a time itself. <3?(Xt) denotes the terminal payoff where <3?(x) is a deterministic 
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function of x. The following approximation procedures can be applied in the same way 
also in the presence of coupon payments. Z and 7 take values in R r and R dxr respectively, 
and " •" in front of the dWt represents the summation for the components of r-dimensional 
Brownian motion. Throughout this paper, we are going to assume that the appropriate 
regularity conditions are satisfied for the necessary treatments. 

2.2 Perturbative Expansion for Non-linear Generator 

In order to solve the pair of (Vt,Zt) in terms of Xt, we extract the linear term from 
the generator / and treat the residual non-linear term as the perturbation to the linear 
FBSDE. We introduce the perturbation parameter e, and then write the equation as 



dV t [e) = c(X t )V} e) dt-eg(X t ,V t (e \zl e) )dt + zl e) ■ dW t (2.4) 

= *(Xr), (2.5) 

where e = 1 corresponds to the original model by0 

f(X t , V t , Z t ) = -c(X t )V t + g{X u V t , Z t ) . (2.6) 



Usually, c(Xt) corresponds to the risk-free interest rate at time t, but it is not a necessary 
condition. One should choose the linear term in such a way that the residual non-linear 
term becomes as small as possible to achieve better convergence. A possible linear term 
9{X)Z in the driver / can be absorbed by the measure change and hence the simple 
reinterpretation of the drift term of the forward components 70 results in the form (|2.4[) . 
See also the discussion in Appendix. 

Now, we are going to expand the solution of BSDE (|2.4p and (|2.5p in terms of e: that 
is, suppose V t {e) and Z[ e) are expanded as 

V^> = + eV« + e 2 vf> + ■ ■ ■ (2.7) 

Z f = Z f) +cZ (i) +£ 2 Z W + .... (2.8) 

Once we obtain the solution up to the certain order, say k for example, then by putting 
e = l, 

V t = j2V®, Z t = Y,zf ] (2.9) 

i=0 i=0 

is expected to provide a reasonable approximation for the original model as long as the 

(i) 

residual term is small enough to allow the perturbative treatment. As we will see, V t 

(i) 

and Z\ , the corrections to each order can be calculated recursively using the results of 
the lower order approximations. 

1 Or, one can consider e — 1 as simply a parameter convenient to count the approximation order. The 
actual quantity that should be small for the approximation is the residual part g. 
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2.3 Recursive Approximation for Perturbed linear FBSDE 
2.3.1 Zero-th Order 

For the zero-th order of e, one can easily see the following equation should be satisfied: 



(0) 



c(X t )V t (0) dt + Zl u> -dW t 



(o) 



(0) 



It can be integrated as 



Vn 



V t {0) = E 



$pf T ) 



e - .17 c ( x °) ds $(x T ] 



(2.10) 
(2.11) 

(2.12) 



which is equivalent to the pricing of a standard European contingent claim. 
Since we have 

e~ So c ( x s) ds $(X T ) = V {0) + / e~ So 4Xs)ds z (o) . m ^ ^ ^ 

J o 



it can be shown that, by applying Malliavin derivative V t , 

j-T 



(e-So< x s) ds <$>(X T )} = V t (e~ SZ < x °) ds ■ dW u + e~ So c{X s )ds z (0) _ 

J t 

Thus, by taking conditional expectation E[- \ Ft], we obtain 



Z<°> = E 



V t (e~ Si c(x B )ds$ 



(Xt)) 



(2.15) 



2.3.2 First Order 



Now, let us consider the process V(0 - y(o). One can see that its dynamics is governed 
by 



d(yW - y/°)) = c(X t )(vf> - - eg(X t , V*>,Z?>)dt + {zf - Z?>) • dW t 
Y^-VP = 0. (2.16) 
Now, by extracting the e-first order terms, we can once again recover the linear FBSDE 



dV t 



(i) 



c(X t )V t m dt - g(X t , V t (U >,ZP)dt + Z?> ■ dW t 



(0) ~(0K 



which leads to 



v±» = o, 



V« = E 



e-s:^ g{Xu ^\z^)du 



Tt 



(2.17) 
(2.18) 

(2.19) 



straightforwardly. By the same arguments in the zero-th order example, we can express 
the volatility term as 

cT 

(2.20) 



ZP = E 



2>t(j[ e-S:^g{X u ,V^,Z^)du 



Tt 



From these results, we can see that the required calculation is nothing more difficult than 
the zero-th order case as long as we have explicit expression for V"(°) and Z^°\ 
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2.3.3 Second and Higher Order Corrections 

We can proceed the same way to the second order correction. By extracting the e-second 
order terms from — (V^ + eV^), one can show that 

dV^ = c(X t )V^dt - (§-g(X t , V<®,ZW) V U + V z g(X t , tf \ Z<°>) • dt + Z« • cW* 



vf> = 



(2.21) 



(2) 

is a relevant FBSDE, which is once again linear in V t . As before, it leads to the following 
expression straightforwardly: 



,(2) 



E 



(2.22) 

T t • 
(2.23) 



In the above calculation, we have assumed the driver function is differentiable. If this is 
not the case, we need to approximate it using some smooth function or apply integration- 
by-parts technique for generalized Wiener functionals (e.g. a composite functional of Dirac 
delta fucntion and a smooth Wiener functional). 

In exactly the same way, one can derive an arbitrarily higher order correction. Due 
to the e in front of the non-linear term g, the system remains to be linear in the every 
order of approximation. However, in order to carry out explicit evaluation, we need to 
give Malliavin derivative explicitly in terms of the forward components. We will discuss 
this issue in the next. 

2.4 Evaluation of Malliavin Derivative 

Firstly, let us introduce a d x d matrix process Y t:U , for u € [t, T], as the solution for the 
following forward SDE: 

d 

d(Y ttV )) = ^(M(^)(^)"^ + ^y(^)(^, u )rrfW- u ) (2.24) 



(Yt 



t.t) 



k=i 
5) , 



(2.25) 



where denotes the differential with respect to the fc-th component of X, and <5*- denotes 
Kronecker delta. 

Now, for Malliavin derivative, we want to express, for u € [t, T], 



E 



V t (e-K< x ^ ds G{X u )) 



(2.26) 
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in terms of Xt, where G is a some deterministic function of X, in general. Thank to the 
known chain rule of Malliavin derivative, we have 



a 

V t (e-K c ^ ds G{X u )) = ^{e-/" c ^^G(I B )(P ( 4) 

>G(X U ) (J U d iC (X s )(V t Xi)ds)y (2.27) 



-J t u c(X s )ds, 



Thus, it is enough for our purpose to evaluate (V t X u ). Since we have 

^ X ^ = jXj t d kl l{X s ){V t X k s )ds + d kl l {X s )(V t X k s )-dW s ) +1 \X t ) (2.28) 

it can be shown that T>tX u follows the next SDE: 

d 

d{V t X l u ) = Y,( d kli( x u){V t Xl)du + d kl \X u ){V t X*)-dW u ) (2.29) 

k=l 

(V t Xl) = f(X t ). (2.30) 

Thus, comparing to Eqs. (|2.24p and (|2.25p . we can conclude that 

(V t X u ) = Y tyUl {X t ) . (2.31) 

As a result, combining the SDE for Yt )U and the Markovian property of X, one can confirm 
that the conditional expectation 

T> t (e-K c{X ° )ds G(X u )) T t ] (2.32) 



E 

is actually given by a some function of Xt. Therefore, in principle, both of the backward 
components can be expressed in terms of Xt in each approximation order. 

In fact, this is an easy task when the underlying process has a known distribution. In 
the next section, we present two such models, and demonstrate how our approximation 
scheme works. We will also compare our approximate solution to the direct numerical 
results obtained from, such as PDE and Monte Carlo simulation. However, in more generic 
situations, we do not know the distribution of X. We will explain how to handle the 
problem in this case using the asymptotic expansion method for the forward components 
in SecJH 



3 Simple Examples 

3.1 A forward agreement with bilateral default risk 

As the first example, we consider a toy model for a forward agreement on a stock with 
bilateral default risk of the contracting parties, the investor (party-1) and its counterparty 
(party-2). The terminal payoff of the contract from the view point of the party-1 is 

<f>(S T ) = S T -K (3.1) 
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where T is the maturity of the contract, and K is a constant. We assume the underlying 
stock follows a simple geometric Brownian motion: 

dS t = rS t dt + aStdWt (3.2) 

where the risk-free interest rate r and the volatility a are assumed to be positive constants. 
The default intensity of party-i hi is specified as 

hi = A, h 2 = X + h (3.3) 

where A and h are also positive constants. In this setup, the pre-default value of the 
contract at time t, Vt, follows H 

dV t = rV t dt- himax(-V t ,Q)dt + h 2 max(y t ,0)dt + Z t dW t 

= (r + X)V t dt + hmax(V t ,0)dt + Z t dW t (3.4) 
V T = $(St) ■ (3.5) 

Now, following the previous arguments, let us introduce the expansion parameter e, 
and consider the following FBSDE: 

dV t (e) = nV t (e) dt-eg(V t {e) )dt + Z { t e) dW t (3.6) 

V T e) = $(St) (3.7) 
dS t = S t (rdt + adW t ) , (3.8) 

where we have defined fi = r + A and g(v) = —hvl{ v>0 y. 

3.1.1 Zero-th order 

In the zero-th order, we have 

dV t {0) = fiV t i0) dt + zl 0) dW t (3.9) 

= §{S T ) • (3.10) 



Hence we simply obtain 



y/ 0) = E e-^ T -V<S>(S T ) 



and 



= e"^-*) (s t e< T ^ - K) (3.11) 
z (0) = e -X(T-t )(jSt (3 12) 



2 See, for example, [3l [5] . 
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3.1.2 First order 

In the first order, we have 

dV 



Thus, we obtain 



V 



(i) 



E 



(i) 
(i) 



fiV t {1) dt - g(V t (0) )dt + Z^dWt 
. 



r 



- f i{T-u) l 



Tt 



= - e -tfr-t) h E\max(S u e r ^~ u ^ -K,0) 



du 



- l) hJ^ C(u;t,S t )du , 



where 



C(u;t,S t ) = S t e r ^N(d 1 {u;t,S t ))-KN(d 2 (u;t,S t )) 
1 / /Se^'K 1 \ 



(3.13) 
(3.14) 



(3.15) 
(3.16) 
(3.17) 

(3.18) 
(3.19) 



and N denotes the cumulative distribution function for the standard normal distribution. 
We can also derive 



iV(di(«;*,S t ))rf« • 



(3.20) 



3.1.3 Second order 



Finally, let us consider the second order value adjustment. In this case, the relevant 
dynamics is given by 



dV 



(2) 



(2) 



^dt-^giy^dt 
. 



As a result, we have 



V 



(2) 



E 



' J\~^ {^ 9 {V^)V^du 



(3.21) 
(3.22) 

(3.23) 



= e-^-^h 2 I I E 

It Ju 



T r T 



1 {S u er(T-u)_ K > 0} C{s;u, S u ) F t dsdu (3.24) 



which can be evaluated as 



pT rT poo 

/ / / 4>(z)C(s;u, S u (z, St))dzdsdu , (3.25) 

J t Ju J —do.(u:t,St) 



't Ju J —d,2(u;t,St) 
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where we have defined 
and 



S u (z,S t ) = £ te (r-^)(«~*)W«=t* 

1 



-±z 2 

e 2 



(3.26) 
(3.27) 



3.1.4 Numerical comparison to PDE 

For this simple model, we can directly evaluate the contract value Vt by numerically solving 
the PDE: 



a 



d 



^-V(t, S) + [rS^-V(t, S) + la 2 S 2 ^V(t, S)) - [ M + hl {v{t>s) > 0} ] V(t, S) = 

with the boundary conditions 

V(T,S) = S-K 

V(t, M) = e-^+W-t) (Me r{T - t] -K), M » K 

V(t, m) = e-^-V (me r(T - 4) -K), m < K . 



PDE and Recursive Approximation 




■8.0 



-10.0 



-12.0 



-PDE 

-1st Order 
-2nd Order 



Maturity (year) 

Figure 1: Numerical Comparison to PDE 

In Fig. [H we have plot the numerical results of the forward contract with bilateral 
default risk with various maturities with the direct solution from the PDE. We have used 



r = 0.02, A = 0.01, h = 0.03, 
a = 0.2, So = 100 , 



(3.28) 
(3.29) 



where the strike K is chosen to make V = for each maturity. We have plot 

for the first order, and + for the second order. Note that we have put e = 1 to 

compare the original model. One can observe how the higher order correction improves the 
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accuracy of approximation. In this example, the counterparty is significantly riskier than 
the investor, and the underlying contract is quite volatile H. Even in this situation, the 
simple approximation to the second order works quite well up to the very long maturity. 

3.2 A self-financing portfolio with differential interest rates 

In this subsection, we consider the valuation of self-financing portfolio under the situation 
where there exists a difference between the lending and borrowing interest rates. Here, we 
consider the problem under the physical measure. 

The dynamics of the self-financing portfolio is governed by [1] 



where r and R are the lending and the borrowing rate, respectively. 9 = (fi — r)/a denotes 
the risk premium. For simplicity, we assume all of the r, R, \i and a are positive constants. 
Here, Zt/a represents the amount invested in the risky asset, i.e. stock St- Let us choose 
the terminal wealth function as 



This spread introduces both of the lending and borrowing activities, which makes the 
problem more interesting. The setup explained here is in fact exactly the same as that of 
adopted by Gobet et al. (2005) [6]. They have carried our detailed numerical studies for 
the above problem and evaluate Vq by regression-based Monte Carlo simulation. In the 
following, we will apply our perturbative approximation scheme to the same problem and 
test its accuracy. 

As usual, let us introduce the expansion parameter as 




(3.30) 



(3.31) 
(3.32) 



$(S T ) 



max(5r — Ki, 0) — 2 max(5T — K% , 0) . 



(3.33) 



rVf'dt - eg(Vf } , Zp)dt + zf dW t 



(3 
(3 



34) 
35) 




where we have defined the non-linear perturbation function as 




(3 



36) 



Now, we are going to expand V t in terms of e. 



3.2.1 Zero-th order 



In the zero-th order, the BSDE reduces to 




(3 
(3 



37) 
38) 



3 Of course, people rarely make such a risky contract to the counterparty in the real market. 
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which allows us to obtain 

V t {0) = E 



{{St -K 1 )+-2{S T -K 2 ) + } 
e-^ T - t HC(K 1 ,S t )-2C(K 2 ,S t )) , 



where we have defined 

C(Ki,S t ) 



d m {Ki,St) = 



S t e^ T - t ^N{d 1 {K l ,S t )) - K t N(d 2 (K t ,S t )) 

)±^(T-t)\ 



In 



for i £ {1,2}. The volatility term is given by 

z f) = e (»-r)(T-t) aS ^ N ( di{KuSt) } _ 2N (d 2 (K 2 ,S t ))) 



(3.39) 
(3.40) 

(3.41) 
(3.42) 

(3.43) 



3.2.2 First order 

Now, in the first order, we have 



dV t 

Vn 



(1) 
(1) 



rV t {1) dt - g(V t {0 \z^)dt + Z^dW t 
. 



(3.44) 
(3.45) 



As before, we can easily integrate it to obtain 

f-T 



V™ = E 



(3.46) 



I e-^-%(V^,Z^)du 
Now, using the zero-th order results, one can show 

g(v} 0) ,Z^) = e-^ T - t \R-r){K l N{d 2 {K l ,S t ))-2K 2 N{d 2 (K 2 ,S t ))} + 

e- r ( T - t \v-r)S t e^ T -^{N(d 1 (K 1 ,S t ))-2N(d 1 (K 2 ,S t ))} (3.47) 

which leads to 



(i) 



e-'P-Q jT dulE 



(R-rj^Nid^K^Su)) - 2K 2 N(d 2 {K 2 ,S u ))) 

- {ix-r)S u e^ T -^(N(d 1 (K 1 ,S u ))-2N{d 1 {K 2 ,S u )))\ Ji] } . 

(3.48) 



By setting 



S U (S U z) = S t e^ u -V exp (~^ 2 (^ - t) + ay/^tz^j 



(3.49) 



we can write the first order correction as 



yW = e ~r(T-t) I du I 

(R - r) [K 1 N{d 2 {K 1 ,S v (S t , z))) - 2K 2 N(d 2 (K 2 , S u (S t , z))) 

(H - r)S u (S t ,z)e^ T - u ^ [n^K^ S u (S t , z))) - 2N(d 1 (K 2 , S u (S t , z)))] } . 

(3.50) 
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The volatility term can also be derived easily as 
d 

n 



I du I dz<j)(z 
Jt Jr 

X (u, z) ^— (KtfickiKi-, S u (S t , z))) - 2K 2 </>(d 2 (K 2 , S u (S t , *)))) 
-( f i-r)S u (S u z)e^ T -^[a( y N(d 1 (K 1 ,S u (S t ,z)))-2N(d l (K 2 ,S u (S t ,z))) 

l={<f>(d 1 (K 1 ,S u (St,z)))-2<l>(d 1 (K 2 ,S u (S t ,z))))}] , (3.51) 



VT-u 
where we have defined 

X (u,z) = 1 if K 1 N{d 2 (K 1 ,S u (S t ,z)))-2K 2 N(d 2 {K 2 ,S u {S t ,z)))>0 

= otherwise . (3.52) 

3.2.3 Second order 

Finally, in the second order, the relevant FBSDE is given by 



dV t 



(2) 



rv^dt- ( g 9 (vr,znvr + f z g(vr,orw j ■>> + zr<™ , 



r(0) ^hT/W 



KO) 7 (0)w(l) 



,(2), 



(2) 



. 



Using the fact that 



we obtain 



9^9(v,z) 

g~ z 9(v,z) 



( T -*) f du f ds f dz! [ dz 2 cp{ Zl )<p{z 2 ) 
Jt Ju Jr Jr 

X (u, zi){-(R - r) 2 [K l N{d 2 (K l ,S s (u, z))) - 2K 2 N(d 2 (K 2 , S s {u, z))) 
+(R -r)(n- r)S s (u, z)e^ T ~^ (jV(di(tfi; S s (u, z))) - 2N(d l (K 2 , S s (u, z)))) } 

x(u, Zl )-ej { 



-(R-r)l { z_ v > 0} 
R-r 

— —hf-vm- 6 



(3.53) 

(3.54) 
(3.55) 



R-r 



(R - r) X (s, u, z)-^= (^^(^(^x, 5 s (n, z))) - 2K 2 <t>{d 2 {K 2 , S s (u, z)))) 
-(ji-^Ssiu^ef^-'HaiNid^KuSsiu,^))) - 2N(d 2 (K 2 , S s (u, z)))\ 



l=[<f>(d 1 (K 1 ,S a (u,^))-2^(d 1 (K 2 ,S 8 (u,z)))]\} 



(3.56) 
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Here, we have defined 

S s (u,z) = S s (S u (S t ,zi),z 2 ) (3.57) 

and also 

X (s,u,z) = 1 if K 1 N{d 2 (K 1 ,S s (u,?)))-2K 2 N{d 2 (K 2 ,S s (u,z)))>0 

= otherwise . (3.58) 

If one needs, it is also straightforward to derive the volatility component. 
3.2.4 Numerical comparison to the result of Gobet et al. 

Gobet et al. (2005) [6] have carried out the detailed numerical study for the above problem 
using the regression-based Monte Carlo simulation. They have used 

^ = 0.05, cr = 0.2, r = 0.01, R = 0.06 

T = 0.25, S = 100, K x = 95, K 2 = 105 . (3.59) 

After trying various sets of basis functions, they have obtained the price as Vq = 2.95 with 
standard deviation 0.01. 

Now, let us provide the results from our perturbative expansion. We have obtained 

y (0) = 2.7863 
y {1) = 0.1814 
V^ 2) = -0.0149 

using the same model inputs. Thus, up to the first order, we have Vq + Vg = 2.968, 
which is already fairly close, and once we include the second order correction, we have 
X^=o = 2.953, which is perfectly consistent with their result of Monte Carlo simula- 
tion. Note that, we have derived analytic formulas with explicit expressions both for the 
contract value and its volatility. 

4 Application of Asymptotic Expansion to Generic Marko- 
vian Forward Processes 

In this section, we consider the situation where the forward components {X{\ consist of 

(i) (i) 

the generic Markovian processes. In this case, we cannot express and Z\ in terms of 
{Xt} exactly, which prohibits us from obtaining the higher order corrections in a simple 
fashion as we have done in the previous section. 

However, notice the fact that what we have to do in each order of expansion is equiva- 
lent to the pricing of generic European contingent claims and hence we can borrow known 
techniques adopted there. In the following, we will explain the use of asymptotic expan- 
sion method, now for the forward components. Although it is impossible to obtain the 
exact result, we can still obtain analytic expression for (V^ , zj: ) up to a certain order 
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of the volatilities of {-Xt}. For the details of asymptotic expansion for volatility, please 
consult with the works QTJ 021 EKE], for example. 

Let us introduce a new expansion parameter 5, which is now for the asymptotic ex- 
pansion for the forward components. We express the relevant SDE of generic Markovian 
process £ M d as 



dX& = l0 (Xi s \S)du + la (xi 5 \5)dW: . 



(4.1) 



Here, we have used Einstein notation which assumes the summation of all the paired 
indexes. For example, in the above equation, the second term means 



la {X<i\8)dW* =J2ia(xi 5 \s)dwz 



(4.2) 



a=l 



We assume 



7a(x,0) = (4.3) 
for a = {1, • • • , r}. Intuitively speaking, it suggests that 5 counts the order of volatility. 



Suppose that, in the (i — l)-th order of e, we succeeded to express V t {i ~ 1] and Z\ 
in terms of xf . Then, in the next order, we can express the backward components as 



r(<-l) 



V, 



z 



(0 



E 
E 



K< x ^ d °G(X( s \6)du 



T 



(4.4) 
(4.5) 



with some function G. If there is no need to obtain (V t 



(i+i) 



, we can just run Monte 



Carlo simulation for X^ to evaluate these quantities in a standard way. However, if we 
want to obtain higher order corrections, we need somehow to express the {y^\ zf^) in 
terms of xf' . 

What we are going to propose here is to expand the backward components around 
5 = 0: 



(0 
(0 



V ^> + sv yj + 5 2 v; l ' 2 > + o(s 2 ) 



Z^ + 5Z^+5 2 Z^ 



+ o{5 2 



(4.6) 
(4.7) 



and express each V t and Z\ in terms of X\ up to a certain order "j" of 5. Although 
we can proceed to arbitrarily higher order of 8, we will present explicit expressions up 
to the second order in this paper. For the interested readers, the work [13] provides the 
systematic methods to obtain higher order corrections. 



Thank to the well-known chain rule for Malliavin derivative, what we have to do is 
only expanding the two fundamental quantities, X u and T>tX u for u £ [t, T], in terms of 
5. Firstly, let us introduce a simpler notation, 



dX® = l0 (X^\5)du + la (X { i\5)dW: 



(4- 
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where a runs through to r with the convention w® = u and = for a G {1, • • • , r}. 

We set the time t- value of as x. Thus our goal is to express and as 
functions of x. We first introduce a d x d matrix process defined as 

d(Y t f)) = d k f a (xi s \6)(Y t ^ d< (4.9) 

(X t f)i = 5], (4.10) 
where d k denotes the differential with respect to the k-ih component of X. Since we have 

l-U l-U 

(X^) i = x i + 1 i > (X^,S)ds+ 1 l(X ( f\5)dW^ (4.11) 
Jt Jt 

applying a Malliavin derivative V t ^ with j3 G {1, • • • , r} gives 

V t>f) (X<Py= / d kl i{xf\8)V t ,p{X^) k ds+ d kl i(xi s \5)V t ,p(xW) k dW« + 7 i(x,5) . 
Jt Jt 

Thus one can show that 

d{V tfi (X^y) = d kl i(xi s \5)V tt e(X^) k dwZ (4.12) 

V^X^y = 7 < (*,*). (4.13) 

Therefore, for a G {1, ■ ■ ■ , r}, we conclude that 



(4.14) 



which implies that the asymptotic expansion of T>tX^ can be obtained from that of 
Therefore, in the following, we first carry out the asymptotic expansion for X and Y. 

4.1 Asymptotic Expansion for and 

We are now going to expand for u G [t, T] as 



and 
where 



X^=xM+SD ttU + ±S 2 E ttU + o(6 2 ) 



Y$=Y ttV + 8H ttV + o(6) , 



Et, u — 



d 2 x u 5) 



and 



5=0 



(0) 



85 2 



5=0 



Yt tU — Y t - U , 



H 



8Y t 



(S) 



t,u 



t,u 



85 



(4.15) 
(4.16) 

(4.17) 
(4.18) 



5=0 
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4.1.1 Zero-th order 

Since 7 (-, 0) = for a 6 {1, • • • , r}, we have 

cLXW = l0 (Xl°\0)du (4.19) 
d(Y t>u )) = d k f (X^ ,0)(F t)tl )}du (4.20) 

with the initial conditions = x and (YtjYj = 5]-, which allows us to express Xl 0) and 
Yt tU as deterministic functions of x. It is also convenient for later calculations to notice 
that Y~ l is the solution of 

d(Y t ^)) = -(Y^id^X^^du (4.21) 

with (Y^)) = 5). 

4.1.2 First order 

By applying d$, we can easily obtain 

d(d s x<f>y = (d^ a (xi s \s)p s (x^ydw^ + d sl i(xi s \d)dw^ (4.22) 
d(d s Y$)i = {47;(4 5 \5)(a,y^ ) ) J fc + a fcn ;(4 5) ^)(a^( 5 ))'(y^) J fc 

+dM(x^ MY^))}dw a u . 

(4.23) 

Putting 5 = 0, they leads to 

dD\ >u = d j7 i(xL \0)Di u du + d 5 f a (xi \0)dwZ (4.24) 

d(H t>u )i = (5 fc7 S(4 \o))(^,jf^ + aH7o(4 \o))z)[ u (y t , u ) J fc ^ 

+ (d k dsf a (xi°\0))(Y t ^dw: . (4.25) 
Now, by using Eq. (|4.2ip . one can show that 

D\, u = (Y t , u )i J\Y t -}){d n k a {s)dw« (4.26) 

(H tjU )) = (Y t<u )l [{Y^fi {(d mn7 i(s))D^ s {Y tjS )f ds + {d m d sl l a {s)){Y t ^dw^) 

(4.27) 

where we have defined the shorthand notation that 

7 j(s):=7j(*i 0) ,0), (4.28) 
which will be used in the following calculations, too. 
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4.1.3 Second order 

Applying dj to the SDE of X® gives us 



d(djx^y 



(a 7a (xW,5))}a|(xf)^% 7 i(xW,5)o a (4 a ))^ (S (^); 

-2d J d sl ^\5)d 8 (X^y +dhl{X^\5)}dwZ . 



Thus, putting 5 = 0, we obtain 

dE\ >u = (5 7 o(x(°),0))j^ u ^ + ^7^(4 \o) J D^A fc u ^ 

+2d J d s fjxi \o)Di u dwZ + d]y a (xi \o)dw: . 

Now we can integrate it as 
El, u = (F t) „)5 f\Y t - 1 ){{d lml k Q (s)DlD^ s ds + fadsi^Dl + d<} 



We do not need the second order terms for 

4.2 Asymptotic Expansion for Malliavin Derivative: V t Xu^ 
For convenience, let us define 

{X*)$ = {V t X^% 

and its expansion as 

id) 1 



where 







5=0 



o2 

/yi\( 2 ) _ _/yi\( j ) 



5=0 



Note that, the zero-th order term (°) vanishes, due to the assumption ()4.3 
From (14,14p . we can easily show that 

mii = (Y ttU ))(d sl i(x,o)) 

mli = (n«)i(^(^o)) + 2(fl t)tl )j(a^(x,o)). 

4.3 Asymptotic Expansion for 

Now, we try to express 

fT 



(i,S) 



E 



e-K< x ^ ds G{X<f,5) 



Ft 



du 
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as a function of x = xf" 1 using the previous results. For that purpose, we first need to 
carry out asymptotic expansion for 

Hg>:= e-K< x( ^ ds G(xi 5 \6) 

to obtain 



<2 = <2+Kn + ^X«+°(<n> 



2o(2) 



where 



<5=0 



(4.38) 
(4.39) 

(4.40) 



<5=0 



Then, we can take the conditional expectation straightforwardly. 



4.3.1 Zero-th order 

We have 



flW = e -jrw)* G(X (o) )0) 



(<>) \ 



(4.41) 



which is a deterministic function of x. 

4.3.2 First order 

It is easy to obtain 



ft u < x ^) ds ^-G{x^ s \s) j u dicix^XdsX^fds 

+(d t G(xi s \5))(dsX^y + dsG(xi 6 \5)} , (4.42) 



which leads to 



^ da {(a i G(X^,0))Di iU + dsG(xi o \0) - G(X^,Q) j U d iC {X^)D{ s ds). 

J t 



(4.43) 



4.3.3 Second order 

In the same way, we can show that 



4 2 1 = e-^ 0) ^{(^G(X^ 

+(^G(xi°\0))El u + d 2 5 G(xi \0) 

+ (T d iC {xf))Dl s d^ G(xi°\0) 

-2([ U d iC (xM)Di s ds^ (^G(X(°),0)4 U + d s G{xV>,0j) 

[d ijC {X^)DlM,s + dA^)E\^ . (4.44) 
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4.3.4 Expression for V t {i ' S) 

Evaluation of the conditional expectation can be easily done by simply applying Ito- 
isometry. Let us first define 



D 



t.u 



(Yt,u)) I {Y t - l ){(dsl k (s))ds 



(4.45) 
(4.46) 

and then we have D\ u = D t u + D\ u . Since the first one is a deterministic function, we 
have, for \/u, s > t, 



E 



DIM 



Dl u Dl s + Dl u Dl, 



where 



uAs 



{Yt, u )i{Y t , s )l I (Y t - l )i(Y^)l(ds7T(v))(ds 7 :(v))dv . 



ui 



(4.47) 
(4.48) 



Then, similarly, we can express 

\Y t ; s y k {di m ^(s)DlbJ s + 20^(8)01 + dh^)}ds . (4.49) 



E t,u 


•■= e i E i 




(Yt,u) 




Using these results, 


R i°l '■- E 


'p(O) 


.Ft 


Rt,l '■= E 


"o(l) 


*] 




and also 






< := E 


"o(2) 


F, 



e -r c ^ 0) )^G(4°),0) 



(4.50) 



(4.51) 



e J t 



J t u c(xi 0) )ds 



{ 



%G(4°),0)^ )UJ D( U + 2(^G(4°>,0))£>; iU + ^G(X(°),0)4 U + 3 5 2 G(4°\0) 
<9 lC (Xj ))^ c (X(°))Z)iX^ 



-25 i G(X(°), ) jf d iC {xf))Di s Diuds 

-2d s G(xi°\0) jy % c(X^)D\ s d S 

-G(xi°\o) J\d ijC {X^)bjM,s + dAX^^ds). 



(4.52) 
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Now, we are able to express V^'^ as a function of x up to the second order of 5 as desired: 

f-T 



(i,S) 



(4.53) 



4.4 Asymptotic Expansion for 



Finally, we are going to express 



r(i,S) 



J t 



=1 E M 



du 



as a function of x = . Let us introduce the two quantities: 



)t,u = e-ft U < x ^ ds d t G(xi 5) ,5)(V t X^y a 

= -e-^ c ^ ) )*G(^,*)(j[ U ^i' ) )(^i* ) )> 

Then, we have 

v t ( e -Jr^>* G( xw *)) = {r] ( a s) ku + ( tf) )t>u . 

Similarly to the previous section, we try to obtain the expressions as 

where both of the zero-th order terms vanish. 
We have 



(4.54) 

(4.55) 
(4.56) 

(4.57) 

(4.58) 
(4.59) 



ds(vi S) )t,u = e~K< x ^ ds {- ^fy iC {xfy){d 5 xfyyd.^ d 3 G{X^\ 5)^X^)1 

+ (d ij G(xi s \s))(d s xWy(v t xWy a + {aid 6 G(x!P,s))(v t xP)i + (d t G{x^\5)){d 5 v t x^y a }. 

(4.60) 



Thus, we obtain 



8=0 



(4.61) 



Similarly we can show that 



8=0 



= -2^fy iC {xfy)Dl^(ri { a\u 

+e - S:<x^) ds [ (ftG(x (0) j ))(^)g + 2(d ij G(^°\ 0))4„(^)g + 2(d i d 5 G(x(°\ 0))(*J)g 

(4.62) 
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In the same way, for ^ , we have 

d s {d)t,u = e-K« x ^ ds G(xW,6) (j\d lC (X^))(d 5 X^yd S 



(d jC (x^))(v t x^y a ds 



- e - R <x?>)d, 



rv 
' Jt 



'{d s G{xjp,6) + (dMx^JWsX^y} (jT d tC (xP)(v t x^y a ds 
c(Jffl)ds G(4^)|" [(d ijC (x^))(d s x^y(v t x^ y a + (d lC (x^))(dsV t xi s) y a 



ds . 
(4.63) 



Thus we can show that 



d 



5=0 



- - e -/" c ^)*G(4 ),0)jr^ftc(xi°)))(^)gcfe (4.64) 



and similarly 



5=0 



-e Jt 



j\d iC {X^))Dlds^ (&\ u 

<^ 0) )*{2 (0,G(X<°>,O) + (^G(X(°),0))AV) j\d jC {xW)){Xif}ds 
G(X<°\0) f [(d^ )))^ . 



(4.65) 



For the evaluation of the conditional expectation, let us define, for Vu > t, 

{H t ,u)) := -E [ ) | ^i] 

= (y titt )» Hy- 1 )* {(d mn7 l (s))Dl s + (d m d sl l (s))}(Y t , s )fds 



(X l a )fl := E 



(i) 



H,u 



Ft 



(Y ttU yAd sl i(x,0)) 



(Y t , u y J (dhl(x,o)) + 2(H tyU y J (dsri(x,o)). 



(4.66) 

(4.67) 
(4.68) 



Using these expressions, one show that 



ft 



e J * 



J t u c(xi 0) )ds 



{(a i G(4 0) ,0))«)S i 1 2-G(4° ) ,0) jf (^°)))(? a )il ) d S } ) 

(4.69) 
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and in the same way that 



(2) , e(2h 



Ft 



= -2 



(^ C (x(°))) J D; iS d S )(f?i i )+d 1) )M 



u « x ? ) rfs {(^G(x(°),o))«)g +2(%G(4°),o)) J D( u (^i)S5 + 2(a,a 5 G(4°),o))«)g 

fcG(4°>,0) + (ftG(4°),0))2r; J f U (d 3 c{xf)))(Xif}ds 



-G(xf,o) jf [(d lC (xi )))(*l)l 2 j + ^^^^ 

Now that we are able to express z[ % ^ as a function of x = x\^ as 

(4^)* = [{mV +t\u + +t\u}du + o(5 2 ) . (4.71) 



(4.70) 



This completes the goal of asymptotic expansion for and Z^ 5 \ which are now 

expressed as functions of x as desired. 

5 Perturbation in PDE Framework 

In this section, we will study the perturbation scheme under the PDE (partial differen- 
tial equation) framework following the so-called four step scheme [8]. We will see that 
our perturbative method makes the four step scheme tractable for the generic situations, 
which only requires standard techniques for the classical parabolic linear PDE. In the next 
section, we will explain the equivalent perturbation method in the probabilistic framework. 

5.1 PDE Formulation based on Four Step Scheme 

Let us consider the following generic coupled non-linear FBSDE: 

dV t = -f(t,X t ,V t ,Z t )dt + ZfdW t 

V T = $(X T ) 

dX t = l0 (t,X u Vt,Z t )dt + ^{t,X t ,V t ,Z t )-dWt 

X = x . (5.1) 

Here, we made the dependence on t explicitly to clearly distinguish it from the stochastic X 
components. As before, we assume that V, Z, X take value in R, IR r and M. d respectively, 
and W denotes a r dimensional standard Brownian motion. 

Following the arguments of the four step scheme of Ma and Yong [9], let us postulate 
that Vt is given by the function of t and Xt as 

V t = v(t,X t ) (5.2) 
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almost surely for Vt G [0, T] . Then, applying Ito's formula, we obtain 

dv t = d t v(t,x t )dt 

+{div(t,X t )4(t,X t ,v(t,Xt),Z t ) + ^d ij v(t,X t )^ i ■ ^){t,X u v(t,X t ),Z t )}dt 
+d iV (t, X t )y% X u v(t, X t ), Z t ) ■ dW t . (5.3) 
Thus, in order that v is the right choice, it should satisfy 
v(T,x) = $(x) (5.4) 
d t v(t,x) + ^div(t,x)io(t,x,v(t,x),z(t,x)) + -d^v^x)^ 1 • 7 J )(t, x, v(t, x), z(t, x))} 

+f(t, x, v(t, x),z(t, x)) = (5.5) 
z(t, x) = diV (t, x)Y(t, x, v (t, x), z(t, x)) , (5-6) 
where the last equation arises to match the volatility term. 

In the four step scheme, one first needs to find the solution z(t, x) satisfying the 
Eq. (|5.6p . And secondly, one has to solve the PDE (|5.5p to obtain v(t, x), which then al- 
lows one to run X as a standalone Markovian process in the third step. And then finally, 
one will obtain the backward components by setting V t = v(t,X t ) and Z t = z(t,X t ). The 
crucial point in the above four step scheme is whether one can finish the step 1 and 2 
successfully. Even if one finds the solution for z, the second step requires to solve the 
non- linear PDE (|5.5p . which is very difficult in general. In the remainder of this section, 
let us study how our perturbation method works to achieve this goal. 



We consider, as before, the original system Eq (|5.ip as a linear decoupled FBSDE with 
perturbations of non-linear generator and feedbacks of the order of e. We write it as 

dV t (e) = c(t, X { t e) )V t (e) dt - eg(t, x£\vf\zfi)dt + Z { t e) ■ dW t 



vp = $p4 e) ) 

dxU = (r(t,xP) + e»(t,xM,vf\zM] 



dt 



+ (a{t, X^ )+e V (t, X® , V® , Z^ )) • dW t 

X^ = x (5.7) 
and the corresponding PDE: 
v^(T,x) = $(x) 

8 t v^ (t, x) + { d iV ^ (t, x) 7 * (t,x, v® , z< e > ) + \dijv^ (t, x) (V • y ) (t, X, , ) } 
+f{t,x,v ( - e \z^) = 

z^(t,x) =d i v (t) (t,x)- f i (t,x,v ( - t \t,x),z^(t,x)) , (5.8) 
where 

f(t,x,v {e \z {e) ) = -c{t,x)v {e) (t,x) + eg{t,x,v {e) {t,x),z {e) (t,x)) (5.9) 

7o(t,x,v (e) ,z (e) ) = r(t,x) + en(t,x,v (e) (t,x),z ( - e) (t,x)) (5.10) 

7(t,x,t; (e) ,z (e) ) = a(t,x) + e??(t, x, u (e) (t, x), (t, x)) . (5.11) 
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We suppose that the solution of the above PDE can be expanded perturbatively in such 
a way that 

v (e) (t,x) = v^(t,x) + ev W (t,x)+e 2 v i2) (t,x) + --- (5.12) 
z^(t,x) = z^(t,x) + ez {1) (t,x) + e 2 z (2 \t,x) + --- , (5.13) 

and then try to solve order by order. If the non-linear terms are small enough, 

we can expect to obtain a good approximation by putting e = 1 in the above expansion 
to a certain order. 

5.2 Zero-th order 

In the zero-th order, the PDE (|5.8|) reduces to 

(d t +C(t,x))v(°\t,x) = 

= $(x) (5.14) 

and 

,x)a l (t,x) . (5.15) 

Here, we have defined the operator C as 

C(t,x) = r i (t,x)d i + -(c/ ■a j )(t,x)d ij - c{t,x) . (5.16) 

This is a standard parabolic PDE and can be handled in the usual way. One can easily 
check that V t = v^(t,X t ) and Z t = z^(t,X t ) solves the FBSDE {577]) when e = 0. 

5.3 First order 

By extracting e-first order terms from the PDE, we obtain 

(d t + C(t, x)) uW (t, x) + (t, x) = 

v^{T,x) = (5.17) 

and 

z (1) (t,x) = d t v^\t,x)a l (t,x) + d i v (0 \t,x)r ) i ^{t,x) . (5.18) 
Here, we have defined 

G {1 \t,x) =d i v i0) (t,x)ij; i(0) {t,x) + d lj v i0) (t,x){a i • V'(°))(t,x) +g {0) {t,x) (5.19) 

and the following notations: 

f^°\t,x) = fi{t,x,v^(t,x),z^(t,x)) (5.20) 
»/ >(t,x) = r)(t,x,v i0) {t,x),z (0 \t,x)) (5.21) 
g<®(t,x) = g(t,x,v^(t,x),z^(t,x)) . (5.22) 

As a result, we once again obtained a linear parabolic PDE. Hence (t)W,z^') can also be 
solved, at least numerically, in a standard fashion. 
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5.4 Second order 

In the second order, one can show v ^ and z^ should satisfy 



(d t + C(t, x))v^ (t,x) + G& (t, x) = 

vW{T,x) = (5.23) 



and 



z^(t,x) = d i v (2) {t,x)a i (t,x) + d i v^{t,x)r ] i ^(t,x) 

+divW (t, x) (vW (t, x)d v + (t, x) ■ V z ) rf^ (t,x) . (5.24) 

Here, is given by 

G^ (t,x) = diV^ (t, x)n 1 ^ (t,x) + div^ (t, x) (v^ (t, x)d v + (t, x) ■ V z ) // , x 

+d lJ v^(t,x)(a i ■ rr>W){t,x) + ^d ijV ^\t,x){rf^ ■ V ^)(t,x) 

+d ij v ( -°\t,x)a i (t,x) ■ (v^(t,x)d v + z^\t,x)-V z )r j j( -°\t,x) 
+(v^{t,x)d v + z^(t,x)-V z )g^(t,x) , (5.25) 

where the partial differentials with respect to v and z are taken by considering fi, r] and g 
as functions of (t,x,v,z). It is still a linear parabolic PDE. 

5.5 Higher orders and an equivalent simpler formulation 

Although we can proceed to higher orders in the same way and solve there is 

another way with a clearer representation. Let us define 

i i 

v®(t,x) = J2 eivU) (t,x)i z®{t,x) = JVztf>(t,z) . (5.26) 

3=0 3=0 

and the operator 

6 k \t,x) = {y (t,x,«W(t,x),^](t,x))^ + i(V-7 m )(t,x,i;W(t,x),zW(t,x))^ ro } 

(5.27) 

Then, one can easily check that the PDE for with (i > 1) can be expressed as 

(d t + rf- 1 } (t, xfj (t, x) + /(*, x, (t, x), z^ (t, x)) = (5.28) 

v®(T,x) = §(x) (5.29) 

and 

z^(t,x) = d lV ^(t,x)j l (t,x,v^(t,x),z^(t,x)) . (5.30) 

It is straightforward to confirm the consistency with the summation of each (v^ k \ z^) for 
(0 < k < i) up to the error terms of o(e*), which is due to the additional e in front of 
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the non-linear terms. Note that, in an arbitrary order, the PDE has a linear parabolic form. 

The above formulation clearly shows that the perturbative treatment of non-linear 
effects of the original system allows us to obtain a series of linear parabolic PDEs with 
the same structure. Solving the PDE for the zero-th order, and then recursively replacing 
the backward components by the solution of the previous expansion order, we can obtain 
an arbitrary higher order of the approximation. 

6 Perturbation in Probabilistic Framework for the Generic 
Coupled Non-linear FBSDEs 

We have now seen the perturbation method in the PDE framework can work even for 
the fully-coupled non-linear FBSDEs. In this section, we will provide a corresponding 
perturbation scheme under the probabilistic framework. As we will see, it is nothing more 
difficult than the decoupled case studied in Sec. O and reduces to the standard calculations 
for the European contingent claims. As a by-product, applying the asymptotic expansion 
method explained in Sec. [U we can also show that it is possible to obtain an analytic 
expression for the non-linear PDE in the Four Step Scheme up to the given order of 
expansion. 

6.1 Generic Formulation 

We try to solve the same FBSDE (15. ip treated in the PDE framework. Suppose that we 
have somehow obtained a solution of (v (t, x), z^~^ (t, x)). Then, let us consider the 
following FBSDE: 

dV® = -f (t,xf ] ,v^(t,xf ] ),z^(t,xf ] ^j dt + zf ■ dW t 

V® = $(a| ] ) 

dxf = l0 (t,xf\v^(t,xP),z^(t,xP))dt 

+7 (t,x¥\vl^t,X®),zW(t,X®)).dW t 
Af = x. (6.1) 

Here, one can immediately check that the solution of the above FBSDE (V^,zf ), as a 
function of (t,xf^), actually satisfies the PDE in the Four Step Scheme given in (|5,28p 
and (|5.30p by setting 

v®(t,x) = V®(t,x), z®(t,x) = Z®(t,x). (6.2) 

Hence the solution of the above FBSDE can be interpreted as the e*-th order approximation 
of the original FBSDE in (|5.ip . Therefore, if we can solve the above FBSDE in probabilistic 
way, we can proceed to an arbitrarily higher order of approximation by simply updating 
the backward components of the non-linear terms recursively. We can also say that it is a 
probabilistic way to solve the non-linear PDE (15. ip order by order of e. 
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One can check that the above FBSDE is actually decoupled and linear by writing it 
explicitly as 

dV® = c(t,xf ] )V t li] dt - eg(t,xf\v^(t,xf ] ),z^(t,xf ] ))dt + zf ■ dW t 

v® = $(x| ] ) 
dxf ] 



r(t, Xf ) + e» (t, Xf , V W (t, X\% Z W (t, xf ))) dt 

+ (a(t,x¥ ] ) + e V (t,xf\v^(t,xM),z^(t,xP))).dW t 



(6.3) 



and hence, we can straightforwardly integrate it as 

r-T 



vr = e 



+ e 



Z li = E 



n 



+ e 



e -St c(s,X*)ds$( X 

e~ J?* < s ' x ^ ds g(u, X® , (u, X® ) , (u, X® ))du 
V t {e-£< s ' x ^ ds ${xf) 

e- /" < s > x ° ] ) ds g(u, XW , (u, X®), z^ (u, X® ))du] 



(6.4) 



(6.5) 



The result is equivalent to the pricing of standard European contingent claims, and 
also has the same form appeared in Sec. [2j Thus, we can apply the asymptotic expansion 
method given in Sec. 2] to the forward components X^ in the same way. This will give us 
the analytical result of (v} l \ zf^) as a function of (t, xf^), up to a given order of volatility 
parameter, say S k . Then we can set 

v®(t,x) = V t [i] (t,x), z®{t,x) = Zf ] (t,x) . (6.6) 
up to the error terms of o(5 k ), and can move on to the higher order of approximations 0. 

6.2 Summary of Recursive Procedures 

Here, let us summarize the procedures of our perturbation method. Firstly, in the zero-th 
order, the corresponding FBSDE is given by 

f ] dt + zf ] -dW t (6.7) 

(6.8) 

t + a{t,xf ] ) -dW t (6.9) 

(6.10) 



dVP 


= <t,X\ 




= S(4 ] 




= r(t,xl 


X [0] 


= x . 



This can be integrated as 



V? 



E 
E 



c( S) 4°]) (is$(x [0] ) 



T 



(6.11) 
(6.12) 



4 Since we finally put 5 = 1 (and also e = 1), the actual order of error terms are of o((a + rf) k ) in this 
example. 



28 



which can be solved either exactly, or analytically up to the certain order of volatility by 
the asymptotic expansion method. Then we set 

«[°](t,x) = V®(t,x) , zM(t,x) = Z®{t,x), (6.13) 

and then put them back in the backward components of (|6.ip with i = 1. We then obtain 
(|6.4p and ()6.5[) with i = 1 . We can express v} 1 ^ and zj^ in terms of t and X^ by using the 
asymptotic expansion method, and use them to define (v^ (t, x), (t, x)) in turn. Now, 
we can move to (|6.4|) and (|6.5|) with i = 2. We repeat the same procedures to the desired 
order of approximation. 

Remark: Although we have considered one-dimensional process for V, it is straight 
forward to extend the method for higher dimensional cases. Once we take the basis of 
X in such a way that the linear drift term V is diagonal, we can proceed without any 
difficulty. The mixing from the other components of V always appears in the lower order 
of e, which keeps the diagonal form of drift term intact in an arbitrary order. 

7 Conclusion and Discussions 

In this paper, we have presented a simple perturbation scheme for non-linear decoupled 
as well as coupled FBSDEs. By considering the interested system as a decoupled linear 
FBSDE with non-linear perturbation terms, we succeeded to provide the analytic approx- 
imation method to an arbitrarily higher order of expansion. We have shown that the 
required calculations in each order are equivalent to those for standard European con- 
tingent claims. We have applied the method to the two simple models and compared 
them with the numerical results directly obtained from the PDE and regression-based 
Monte Carlo simulation. Both of the examples clearly demonstrated the strength of our 
method. We have also shown that the use of the asymptotic expansion method for forward 
components allows us to proceed to the higher order of perturbation even if the forward 
components do not have known distributions. 

In the last part of the paper, we have studied the perturbative method in the PDE 
framework based on the so-called Four Step Scheme. We have shown that our perturba- 
tive treatment renders the original non-linear FBSDE into the series of linear parabolic 
PDEs, which are straightforward to handle. Furthermore, by the equivalence of the two 
approaches, we were also able to provide the corresponding perturbative method in prob- 
abilistic framework which is explicitly consistent with the Four Step Scheme up to a given 
order of expansion. 

The perturbation theory presented in this paper may turn out to be crucial to inves- 
tigate various interesting problems, such as those given in the introduction, which have 
been preventing analytical treatment so far. The application of the new method to the 
important financial problems is one of our ongoing research topics. 

Finally, let us remark on the further extension to the cases including jumps. Although, 
in this work, we have only considered the dynamics driven by Brownian motions, the 
same approximation scheme can also be applied to more generic cases. Although it will 
be more difficult to obtain explicit expressions in terms of forward components, if we 
choose the specific forward processes with appropriate analytical properties, we should 
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be able to proceed in the similar way. Particularly, the separation of the original system 
into the decoupled linear FBSDE and the non-linear perturbation terms can be done in a 
completely parallel fashion. 



A Linear Volatility Term in the Driver 

As we have briefly mentioned in Section [21 there are the situations where the driver contain 
the linear term in volatility Although it is possible to absorb it by the measure change, 
there may be the situation where its direct treatment improves the numerical performance 
than inducing the drift modification in the forward components. Consider the situation 
where the backward component follows the following BSDE: 



dV t = c(t, X t )V t + 8(t, X t ) -ZAdt- eg(t, X t ,V t , Z t )dt + Z t ■ dW t . 



(A.l) 



Also in this case, the recursive procedures can be carried out in the same way by simply 
replacing the factor 



cxp 



c(u, X u )du 



by 



£(t,s) =exp (-J (c(u,X u ) + ^\\0(u,X u 



W du 



(u,X u )-dW u 



For example, the results of Section 12.3.11 can now be expressed as 



(0) 



Z 



(o) 



E 
E 



£{t,T)*(X T ) 



V t (£(t,T)$(X T ) 



Ft 



(A.2) 

(A.3) 

(A.4) 
(A.5) 



Higher order and coupled cases can be expressed similarly. In general, it is difficult to say 
which method performs better, since both of them have errors of the same order of e at a 
given expansion. 

B Another expansion method for coupled FBSDEs 

In the appendix, we provide another method which is more closely related to that of SecJ2] 
for generic coupled FBSDEs. We consider the same coupled FBSDE as in Sec. El 



= c(X^)V t {i) dt - eg(Xr,Vr,ZDdt + Zf> • dW t 



r(e) 



(B.l) 
(B.2) 



dX^ = r{Xf>)dt + en{X?>, Vf>, Zf>)di + ( *{X?>) + er,(X?> , Vf>, Zf>) ) • dW t (B.3) 



.(e)^ 



-(e) T/ (e) 7 (e). 



4 £) = - 



(B.4) 



Here, we absorbed a possible explicit time dependency to X. 



J We are grateful for an anonymous referee to point this out. 
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As in Secj2j we are going to expand the backward as well as forward components in 
terms of e. Suppose that we have 



t 



V t {0) + eV t v± > + e z V t ^> + 
Z (0) + eZ W +e 2 z (2) + 



r(l) 



(2) 



X t {0) + eX^ + + • • • 



(2) 



(B.5) 
(B.6) 
(B.7) 



Now, let us derive each term separately. 



B.l Zero-th order 

In the zero-th order, the original equation reduces to a linear and decoupled FBSDE, 
which will going to serve as the center of expansion: 



dV t 



(0) 



c(xP)v} 0) dt + zP.dW t 



V, 



(o) 



T 



$(X 



(oh 

T i 



X 



(o) 
t 

(0) 



r(X t (0) )cft + ^(X-, 



(Oh 



dWt 



x . 



(B.8) 
(B.9) 
(B.10) 
(B.ll) 



x(°) is now the standard Markovian process completely decoupled from the backward 
components. We can easily solve the backward components as 



(0) 



(0) 



E 
E 



e - It c ( x s 0) ) ds $(X, 



V t (e-f'< x ^ d ^(X^) 



(B.12) 
(B.13) 



As explained in the previous section, we can express V^ and zf^ in terms of xf^ by the 
help of the asymptotic expansion for volatility even if the process of X^ ) does not have 
known distribution. 



B.2 First order correction 

Now let us consider the dynamics of V^ — V^ and X^ — X(°> . Following the same 
arguments in Sec 121 one can show straightforwardly that 



dV, 



(i) 



= c(xi 0) )V t (1) dt - g w {t)dt + Z\ Li ■ dW t 

$W(r) 



(i) 



v. 



(i) 



T 



dX. 



(i) 



V x r(X t (0) ) • X^'dt + fiW(t)dt + r]W(t) ■ dW t 



(0), 



X 



(i> 







(B.14) 
(B.15) 
(B.16) 
(B.17) 
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(0) T/ (0) z (p)j 



a(xr,v t , 



where we have used shorthand notations: 
$M(T) 

A* (0) (*) 
„<°>(t) 



(B.18) 
(B.19) 



(0)> 



*( i )(4 > > 4 1 >) = v x #(4 0) )-4 1) 

5 (°)(t)-(v, C (4)).4y/°) (B.20) 

(B.21) 
(B.22) 
(B.23) 



m4V/°U (0) ) 

^(0)^(0)^(0)) 

^) ( 4^4V t (0 U (0) 



(i) 



^ )(t) + v x a(4 >; 

(0) 



Note that, since we have obtained V t and Z\ as the functions of X\ , the pair 
consists of a Markovian process, which is indeed decoupled from (V^ 1 ), Z^). 
This means that we have ended up with the decoupled linear FBSDE also for the first 
order correction. Hence, one can easily solve the backward components as 



e -/ t T ^ 0) )^$(i)(T) + J e-K< x ^ ds gW{u)du 
Z W = E ~ Vt ( e -f?c(xW) dS Q(i) {T) + j\-i:<x^) d s g {i) {u)d , 



(B.24) 
(B.25) 



As a result, we have obtained V^as well as as the functions of xj. ^ and x\ l \ 
B.3 Second and Higher order corrections 

We can continue to the higher order corrections in the same way. By considering — 
(y(0) + e y(i)) an d _ + e X^), and then extracting the e-second order terms, 

one can show that 



dV t {2) 

,(2) . 
' T ■ 



c(xl 0) )V t (2) dt - g (2) {t)dt + Z {2) ■ dW t 



dx[ 2) = V x r{Xl u> ) ■ Xf'dt + n^{t)dt + r]W(t) • dW t 



-(0)^. Y Wj + ^ .,m mdt , _(2), 



4 2) = . 



Here we have defined 



$(2)( T ) = $W (4°) ,4), 4 2 ); 



= v x <s>(xW).xp + fa*{xW)xF>xp 
9 v {t) = ^ {x (o) x (d > x m F (o f (d z (o } 

= {v^°)(t) • xf) + + v^°)(t) • z« 



'(2) 



-(O)^yi(l) v j(l) 



(B.26) 
(B.27) 
(B.28) 
(B.29) 

(B.30) 
(B.31) 
(B.32) 



(v x c(X 



, , • 4> + ^ciX^X^xi^V^ - (V.c(4 0) ) • 4 X) ) ^ 



(0), 



(i) 



(B.33) 
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And also for the process , we have defined 



^Xx?\x?\x?\v?\v?\z?\zV) 



and 



+ ^d ij r{xf)X^X^ 



(B.34) 
(B.35) 

(B.36) 
(B.37) 



+ 



{v^°)(i)X t (1) + ^V {0) W t (1) + V zV M(t) ■ Z^} . (B.38) 



One can check that the pair (X t (0) , X[ 1] , Xf ] ) consists of a Markovian process and it 
is decoupled from the backward components. Hence, once again, we have obtained the 
decoupled linear FBSDE, which is solvable as before: 



V, 



(2) 



z 



(2) 



E 
E 



e -/-c(xi°>) ds$ (2) (T) + J e -i:<X^) dSg{ 2) {u)du 



(B.39) 
(B.40) 



In completely the same way, we can proceed to an arbitrarily higher order correction. 
In each expansion order e l , the set {xf 1 , • • • , X}: ) follows Markovian process decoupled 
from the backward components, and also the FBSDE continues to be linear thank to the 
e in front of the non- linear terms. 



Remark: For simplicity of the presentation, we have used the common perturbation param- 
eter e both for the non-linearity in backward components as well as the feedback effects in 
the forward components. However, as one can easily expect, it is also possible to introduce 
multiple perturbation parameters. 



B.4 Consistency to the result of Sec. [6] 

For completeness, let us check the consistency to the result of Sec. [H In the zero-th order, 
the corresponding FBSDEs are exactly equal. In the first order, we have obtained 



V, 



[i] 



E 



+e 



T 



* <'^g (u, X W , v (0) (u , X W ) , z<°> (u, X W )) du 



Ft 



(B.41) 



and z] 1 ' as its Malliavin derivative in Sec. [6l Now, in order to compare it to the current 
method, let us consider 



[i] 



v t 



(0) 



(B.42) 
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and extract e first order terms by expanding xf^ = xf^ + eX^ . Since the second term 
of (|B.41|) has already e in the front, we simply get a contribution from there as 



E 



7i 



(B.43) 



From the first term, we have to expand X™ in the terminal payoff and also in the discount: 



where we have used the Gateaux derivative H. Second term can be interpreted that there 
is a change in value by 

-(V x c(u,X^).X^)V^ (B.44) 

at each point of time u, which is summed and discounted back to the current time. By 
summing the above two terms, one can easily confirm its equivalence to the result of 
the previous section. Applying Malliavin derivative automatically tells the consistency of 
volatility terms. Using the same arguments, one can check the consistency between the 
current method and that of Sec. [H Although we have solved V® separately, the sum 
J2k=o e k V^ can be shown equivalent to V"M up to the error terms o(e l ). 
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